
#include "weightsScheme.h"

inline Number omega(Number d) {
   
   const Number NINFINITY = std::numeric_limits<Number>::infinity();
   
   if (d==NINFINITY) return 0.00000001;
   if (d==0)         return 1e8;
   else 
      return 1.0/(pow(d,2) + 0.00000001);
}

Number* compute_WeightsTable(Mesh *model, const Constraints *constraints) {
   
   cout<<" compute_WeightsTable(Shape3D *model)"<<endl;
   
   uint texSize = (uint)(std::ceil(std::sqrt(model->nvertices())));
   uint texSizetexSize = texSize*texSize;
   Number *weightsTable = new Number[texSizetexSize*constraints->size()];
   
   for (unsigned j=0; j<constraints->size(); ++j) {
      model->floodfill_Dijkstra(constraints->ids[j]);
      Vertex_iterator vit = model->vertices_begin();
      for (unsigned i=0; i<model->nvertices(); ++i, ++vit) {
         weightsTable[j*texSizetexSize + i] = omega(vit->d);
      }
   }
   
   return weightsTable;
}

Number* compute_WeightsTable_Trihedrons(Mesh *model, const Constraints *constraints) {
   
   cout<<" compute_WeightsTable_Trihedrons(Mesh *model)"<<endl;
   
   uint texSize = (uint)(std::ceil(std::sqrt(model->nvertices())));
   uint texSizetexSize = texSize*texSize;
   Number *weightsTable = new Number[texSizetexSize*constraints->size()];
   
   const Number *a = constraints->getPositionPointer(0);
   const Number *b = constraints->getPositionPointer(1);
   const Vector3 v(Point3(a[0],a[1],a[2]), Point3(b[0],b[1],b[2]));
   const Number d = sqrt(v*v);
   
   for (unsigned j=0; j<constraints->size(); j+=4) {
      model->floodfill_Dijkstra(constraints->ids[j]);
      Vertex_iterator vit = model->vertices_begin();
      for (unsigned i=0; i<model->nvertices(); ++i, ++vit) {
         Number pd = vit->d;
         weightsTable[j*texSizetexSize + i] = omega(pd);
         weightsTable[(j+1)*texSizetexSize + i] = 
         weightsTable[(j+2)*texSizetexSize + i] = 
         weightsTable[(j+3)*texSizetexSize + i] = omega(pd+d);
      }
   }
   
   return weightsTable;
}
